P1829 [国家集训队]Crash的数字表格

$\sum_{i=1}^n\sum_{j=1}^mlcm(i,j)$

$\sum_{i=1}^n\sum_{j=1}^m\frac{i*j}{gcd(i,j)}$

$\sum_{d=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m [gcd(i,j)=d] * \frac{i*j}{d}$

$\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}[gcd(i,j)=1]*i*j*d$

$\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}\sum_{k|gcd(i,j)}\mu(k)*i*j*d$

$\sum_{d=1}^{min(n,m)}\sum_{k=1}^{min(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)}\mu(k)*k^2*d\sum_{i=1}^{\lfloor \frac{n}{dk} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{dk} \rfloor}i*j$

$\sum_{d=1}^{min(n,m)}d\sum_{k=1}^{min(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)}\mu(k)*k^2*(\frac{\lfloor \frac{n}{dk} \rfloor*(\lfloor \frac{n}{dk} \rfloor+1)*\lfloor \frac{m}{dk} \rfloor*(\lfloor \frac{m}{dk} \rfloor+1)}{4})$

内外进行两次数论分块，可以做到$\Theta(\sqrt n * \sqrt n~)=\Theta(n)$求解

#include <cstdio>
#include <iostream>
using namespace std;
#define Mod 20101009

const int MAXN = 10000000;
int n , m , f[ MAXN + 5 ];

int k , prime[ MAXN + 5 ] , mu[ MAXN + 5 ];
bool vis[ MAXN + 5 ];
void sieve( ) {
    mu[ 1 ] = 1;
	for( int i = 2 ; i <= MAXN ; i ++ ) {
		if( !vis[ i ] ) {
			prime[ ++ k ] = i;
			mu[ i ] = -1;
		}
		for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
			vis[ i * prime[ j ] ] = 1;
			if( i % prime[ j ] == 0 ) break;
            mu[ i * prime[ j ] ] = -mu[ i ];
		}
	}
    for( int i = 1 ; i <= MAXN ; i ++ )
        f[ i ] = ( f[ i - 1 ] + 1ll * mu[ i ] % Mod * i % Mod * i % Mod ) % Mod;
}

int Quick_pow( int x , int po ) {
    int Ans = 1;
    while( po ) {
        if( po % 2 )
            Ans = 1ll * Ans * x % Mod;
        x = 1ll * x * x % Mod;
        po /= 2;
    }
    return Ans;
}
int inv( int x ) {
    return Quick_pow( x , Mod - 2 );
}

int Get( int n1 , int m1 ) {
    int d1 = min( n1 , m1 ) , Ans = 0;
    for( int l = 1 , r ; l <= d1 ; l = r + 1 ) {
        r = min( n1 / ( n1 / l ) , m1 / ( m1 / l ) );
        Ans = ( Ans + 1ll * ( f[ r ] - f[ l - 1 ] + Mod ) % Mod * ( 1ll * ( n1 / l ) * ( n1 / l + 1 ) % Mod * ( m1 / l ) % Mod * ( m1 / l + 1 ) % Mod ) % Mod * inv( 4 ) % Mod ) % Mod;
    }
    return Ans;
}
int solve( int n , int m ) {
    int Ans = 0 , d = min( n , m );
    for( int l = 1 , r ; l <= d ; l = r + 1 ) {
        r = min( n / ( n / l ) , m / ( m / l ) );
        Ans = ( Ans + 1ll * Get( n / l , m / l ) * ( l + r ) * ( r - l + 1 ) / 2 ) % Mod;
    }
    return Ans;
}

int main( ) {
    sieve( );
    scanf("%d %d",&n,&m);
    printf("%d\n",solve( n , m ));
    return 0;
}